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Abstract 

The dipole moment of any finite and neutral system, having a square- 
integrable wavefunction, is a well defined quantity. The same quantity 
is ill-defined for an extended system, whose wavefunction invariably 
obeys periodic (Born-von Karman) boundary conditions. Despite this 
fact, macroscopic polarization is a theoretically accessible quantity, for 
either uncorrelated or correlated many-electron systems: in both cases, 
polarization is a rather "exotic" observable. For an uncorrelated — either 
Hartree-Fock or Kohn-Sham — crystalline solid, polarization has been 
expressed and computed as a Berry phase of the Bloch orbitals (since 
1993). The case of a correlated and/or disordered system received a 
definitive solution only very recently (1998): this latest development 
allows us to present here the whole theory from a novel, and very general, 
viewpoint. The modern theory of polarization is even relevant to the 
foundations of density functional theory in extended systems. 



1 Introduction 

The dipole moment of any finite A^-electron system in its ground state is a 
simple and well defined quantity. Given the many-body wavefunction \1/ and 
the corresponding single-particle density n{r) the electronic contribution to 
the dipole is: 

(R) = / dr rn(r) = (^|R|^), (1) 



where R = J2iLi^i (atomic Hartree units are adopted throughout). This 
looks very trivial, but we are exploiting here an essential fact: the ground 
wavefunction of any finite A^-electron system is square-integrable and vanishes 
exponentially at infinity; the density vanishes exponentially as well. 
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Considering now a macroscopic solid, the related quantity is macroscopic 
polarization, which is a very essential concept in any phenomenological 
description of dielectric media [0: this quantity is ideally defined as the dipole 
of a macroscopic sample, divided by its volume. The point is that, when using 
Eq. (HI), the integral is dominated by what happens at the surface of the sample: 
knowledge of the electronic distribution in the bulk region is not enough to 
unambiguosly determine the dipole. This looks like a paradox, since in the 
thermodynamic limit macroscopic polarization must be an intensive quantity, 
insensitive to surface effects. 

Macroscopic polarization in the bulk region of the solid must be determined 
by what "happens" in the bulk as well. This is the case if one assumes a 
model of discrete and well separated dipoles, a la Clausius-Mossotti: but real 
dielectrics are very much different from such an extreme model. The valence 
electronic distribution is continuous, and often very delocalized (particularly in 
covalent dielectrics). Most textbooks attempt at explaining the polarization of 
a periodic crystal via the dipole moment of a unit cell, or something of the kind 
[0, ^. These definitions are incorrect 0]: according to the modern viewpoint, 
bulk macroscopic polarization is a physical observable completely independent 
from the periodic charge distribution of the polarized crystalline dielectric. 

In condensed matter physics the standard way for getting rid of undesired 
surface effects is to adopt periodic Born-von Karman boundary conditions 
(BvK). Indeed, the BvK choice is mandatory in order to introduce even the 
most elementary topics, such as the free-electron gas and its Fermi energy, 
or the Bloch theorem 0, Unfortunately, the adoption of BvK does not 
solve the polarization problem. In fact the dipole cannot be evaluated as in 
Eq. (0) when the wavefunction obeys BvK: the integrals are ill-defined due to 
the unbounded nature of the quantum-mechanical position operator. 

For this reason macroscopic polarization remained a major challenge in 
electronic structure theory for many years. The breakthrough came in 1992, 
when polarization was defined in terms of the wavef unctions, not of the charge. 
This definition has an unambiguous thermodynamic limit, such that BvK 
and Bloch states can be used with no harm 0. In the following months a 
modern theory of macroscopic polarization in crystalline dielectrics has been 
completely established 0, thanks to a major advance due to R.D. King- 
Smith and D. Vanderbilt |p, who expressed polarization in terms of a Berry 
phase H, [ill. A comprehensive account of the modern theory exists 0. 
Other less technical presentations are available as well [p!2| , p!3| , p!^ ; for an 
oversimplified nontechnical outline see Ref. [l^. First-principle calculations 
based on this theory have been performed for several crystalline materials, 
either within DFT (using various basis sets) or within HF (using LCAO basis 
sets) [m. 



2 



All of the above quoted work refers to a crystalline system within an 
independent-electron formulation: the single-particle orbitals have then the 
Bloch form, and the macroscopic polarization is evaluated as a Berry phase 
of them. The related, but substantially different, problem of macroscopic 
polarization in a correlated many-electron system was first solved by Ortiz and 
Martin in 1994 [jl^. However, according to them, polarization is defined — and 



computed 1T|, H]— by means of a peculiar "ensemble average", integrating 
over a set of different electronic ground states: this was much later (1998) 



shown to be unnecessary. In Ref. in fact, a simpler viewpoint is taken: 



the polarization of a correlated solid is defined by means of a "pure state" 
expectation value, although a rather exotic kind of one. By the same token, 
it was also possible to define [0, ^ |2T| — and to compute — macroscopic 
polarization in noncrystalline systems. 

In the present work we take advantage of the most recent developments for 
reconsidering the whole theory under a new light. At variance with previous 
presentations we are not going to introduce explicitly the Berry phase concept: 
in fact, within the formulation of Ref. the Berry phase appears very much 
"in disguise". Instead, a major role is played by a precursor work [|^, 
apparently unrelated either to the polarization problem or to a Berry phase, 
which will be reexamined here and used to introduce the polarization theory. 

Finally, let me just mention a latest development, not to be discussed 
in the present work, where ideas spawned from the polarization theory — 
and more specifically from Ref. — are used to investigate wavefunction 
localization ||25|| . 



2 The "electron in broth" 

Adopting a given choice for the boundary conditions is tantamount to defining 
the Hilbert space where our solutions of Schrodinger's equation live. For the 
sake of simplicity, I am presenting the basic concept by means of the one- 
dimensional case. For a single-particle wavefunction BvK reads ip{x + L) = 
ipi^x), where L is the imposed periodicity, chosen to be large with respect to 
atomic dimensions. Notice that lattice periodicity is not assumed, and BvK 
applies to disordered systems as well. 

By definition, an operator maps any vector of the given Hilbert space 
into another vector belonging to the same space: the multiplicative position 
operator x is therefore not a legitimate operator when BvK are adopted for the 
state vectors, since xil){x) is not a periodic function whenever iIj{x) is such. 
It is then obvious why Eq. (|l|) cannot be used in condensed matter theory. 
Of course, any periodic function of x is a legitimate multiplicative operator in 
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the Hilbert space: this is the case e.g. of the nuclear potential acting on the 
electrons. 

Before switching to the polarization problem, it is expedient to discuss an 
important precursor work, apparently unrelated to the polarization problem, 
where nonetheless the expectation value of the position operator plays the key 
role. 

Some years ago, A. Selloni et al. |^ addressed the properties of electrons 
dissolved in molten salts at high dilution, in a paper which at the time was 
commonly nicknamed the "electron in broth". The physical problem was 
studied by means of a mixed quantum-classical simulation, where a lone 
electron was adiabatically moving in a molten salt (the "broth") at finite 
temperature. The simulation cell contained 32 cations, 31 anions, and a single 
electron. KCl was the original case study, which therefore addressed the liquid 
state analogue of an F center; other systems were studied afterwards . 
The motion of the ions was assumed as completely classical, and the Newton 
equations of motion were integrated by means of standard molecular dynamics 
(MD) techniques, though the ionic motion was coupled to the quantum degree 
of freedom of the electron. The electronic ground wavefunction was determined 
solving the time-dependent Schrodinger's equation at each MD time step. As 
usual in MD simulations, periodic boundary conditions were adopted for the 
classical ionic motion. Ideally, the ionic motion occurs in a simulation cell 
which is surrounded by periodic replicas: inter-cell interactions are accounted 
for, thus avoiding surface effects. Analogously, the electronic wavefunction is 
chosen in the work of Selloni et al. to obey BvK over the simulation cell, 
and therefore features periodic replicas as well. A plot of such an electronic 
distribution, in a schematic one-dimensional analogue, is given in Fig. 1. 



One of the main properties investigated in Ref. |£3| was the electronic 
diffusion, where the thermal ionic motion is the driving agent (within the 
adiabatic approximation). In order to perform this study, one has to identify 
first of all where the "center" of the electronic distribution is. Intuitively, the 
distribution in Fig. 1 appears to have a "center" , which however is defined only 
modulo the replica periodicity, and furthermore cannot be evaluated simply 
as in Eq. (P, i.e. {x) = J dx x\iIj{x)\^, precisely because of BvK. Selloni et 
al. solved the problem by means of a very elegant and far-reaching formula. 



presented below. The work of Ref. |^ can be regarded as the many-body 
generalization of it. 
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Figure 1: The distribution of a single-particle orbital within periodic 

Born-von-Karman boundary conditions 

3 The main formula: One electron 



According to Refs. |20|, the key quantity for dealing with the 



position operator within BvK is the dimensionless complex number 3, defined 
as: 

3 = (V^le^^l^) = tdxe^'-^^mx)]', (2) 
^0 

whose modulus is no larger than 1. The most general electron density, such 
as the one depicted in Fig. 1, can always be written as a superposition of a 
function nioc(x), normalized over {—00, 00), and of its periodic replicas: 

00 

= nioc{x - Xo - mL). (3) 



Both Xo and nioc(x) have a large arbitrariness: we restrict it a little bit 
by imposing that Xq is the center of the distribution, in the sense that 
dxxnioc{x) = 0. Using Eq. (||), 3 can be expressed in terms of the Fourier 
transform of nioc as: 

3 = e^^-ni..(-^). (4) 

If the electron is localized in a region of space much smaller than L, its Fourier 
transform is smooth over reciprocal distances of the order of and can be 
expanded as: 



27r. 1 /27r' ^ 



00 



^^ioc(--^) = 1- 2 ( ^ ) / dxx'n,oc{x) + 0{L-'). (5) 
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A very natural definition of tlie center of a localized periodic distribution 
|?/'(a;)p is therefore provided by the phase of 3 as: 

(a;) = ^Im ln3, (6) 



which is in fact the formula first proposed by Selloni et al. [23, 24 1. The 
expectation value (x) is defined modulo L, as expected since |'?/'(a;)p is BvK 
periodic. 

The above expressions imply {x) ~ Xq mod L; in the special case where 
^ioc(3;) can be taken as an even (centrosymmetric) function, its Fourier 
transform is real and Eq. (^) yields indeed (x) = Xq mod L. In the case 
of extreme delocalization we have instead = 1/L and 3 = 0: hence the 

center of the distribution (x), according to Eq. (|^), is ill-defined. For a more 



general delocalized state, we expect that 3 goes to zero at large L [p5|] . 

We have therefore arrived at a definition of (x) within BvK which has 
many of the desirable features we were looking for: nonetheless, there is a 
property that is even more important, and which we are going to demonstrate 
now. Suppose the potential which the electron moves in has a slow time 



dependence — as was the case in Ref. [23| — and we wish to follow the adiabatic 



evolution of the electronic state If we call \ipj) the instantaneous 

eigenstates at time t, the lowest order adiabatic evolution of the ground-state 



density matrix is |26]: 



where the phases have been chosen in order to make ly^o) orthogonal to its 
time derivative l^o)- The macroscopic electrical current flowing through the 
system at time t is therefore: 

0) = -i-(^W>^-i^ <^°IPl^'><^'l^°' +cc. (8) 

It is then rather straightforward to prove that (j) to lowest order in 1/L equals 
— {l/L)d{x)/dt, where (x) is evaluated using in Eq. (^) the instantaneous 
ground eigenstate: 

~ "^^"^ l^"" (<^o|e'^"|^o)- (9) 

This finding proves the value of the "electron-in-broth" formula, Eqs. (0) and 
(P) in studying electron transport p3, 12^ . 
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4 The main formula: Many electrons 



So much about the one-electron problem: we are now going to consider a finite 
density of electrons in the periodic box. To start with, irrelevant spin variables 
will be neglected: for the sake of notation simplicity, I will first illustrate the 
main concepts on a system of "spinless electrons". 

Even for a system of independent electrons, our approach takes a simple 
and compact form if a many-body formulation is adopted. BvK then imposes 
periodicity in each electronic variable separately: 

"^oiXl, Xi, Xat) = ^0(3^1, ••• , Xi + L, . . . , Xat). (10) 

Our interest is indeed in studying a bulk system: electrons in a segment 
of length L, where eventually the thermodynamic limit is taken: L 00, 
N — * 00, and N/L = hq constant. We also assume the ground state 
nondegenerate, and we deal with insulating systems only: this means that 
the gap between the ground eigenvalue and the excited ones remains finite for 
L 00. 

We start defining the one-dimensional analogue of R, namely, the multi- 
plicative operator X = YliLi complex number 

3^ = (vi/|e^^^|vl/). (11) 

It is obvious that the operator X is ill-defined in our Hilbert space, while its 
complex exponential appearing in Eq. (|11|) is well defined. The main result of 



Ref. 1 20] is that the ground-state expectation value of the position operator is 



given by the analogue of Eq. (^, namely: 

(X) = ^Im In (12) 

a quantity defined modulo L as above. 

The right-hand side of Eq. (|12D is not simply the expectation value of 
an operator: the given form, as the imaginary part of a logarithm, is indeed 
essential. Furthermore, its main ingredient is the expectation value of the 
multiplicative operator e*"^^"^: it is important to realize that this is a genuine 
many-body operator. In general, one defines an operator to be one-body 
whenever it is the sum of N identical operators, acting on each electronic 
coordinate separately: for instance, the X operator is such. In order to 
express the expectation value of a one-body operator the full many-body 
wavefunction is not needed: knowledge of the one-body reduced density matrix 
p is enough: I stress that, instead, the expectation value of e*"?"^ over a 
correlated wavefunction cannot be expressed in terms of p, and knowledge 
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of the A^-electron wavefunction is explicitly needed. In the special case of 
a single-determinant, the A^-particle wavefunction is uniquely determined by 
the one-body reduced density matrix p (which is the projector over the set 
of the occupied single-particle orbitals): therefore the expectation value (X), 
Eq. (|l2l), is uniquely determined by p. But this is peculiar to uncorrected 
wavefunctions only: this case is discussed in detail below. 

As in the one-body case, whenever the many-body Hamiltonian is slowly 
varying in time, the macroscopic electrical current flowing through the system 
is given by 

- -i^(A->, (13) 

where (X) is evaluated using in Eq. (|TTp the instantaneous ground eigenstate 
of the Hamiltonian at time t: this result is proved in Ref. [Q. Considering 



now the limit of a large system, {X) is an extensive quantity: the macroscopic 
current (J), Eq. ([T3|), goes therefore to a well defined thermodynamic limit. 

We stress that nowhere in our presentation have we assumed crystalline 
periodicity. Therefore our definition of (X) is very general: it applies to any 
condensed system, either ordered or disordered, either independent-electron 
or correlated. 



5 Macroscopic polarization 

In the Introduction, we have discussed what polarization is not, by outlining 
some incorrect definitions @, |^ and their problems [Q. We have not stated 
yet what polarization really is: to this aim, a few experimental facts are 
worth recalling. The absolute polarization of a crystal in a given state has 
never been measured as a bulk property, independent of sample termination. 




Figure 2: A meaurement of the piezo- 
electric effect. The current flowing 
along a crystal shorted out by a ca- 
pacitor is measured while the crystal is 
strained in one direction. 
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Instead, well known bulk properties are derivatives of the polarization with 
respect to suitable perturbations: permittivity, pyroelectricity, piezoelectricity, 
dynamical charges. In one important case — namely, ferroelectricity — the 
relevant bulk property is inferred from the measurement of a finite difference 
(polarization reversal). In all cases, the derivative or the difference in 
the polarization is typically accessed via the measurement of a macroscopic 
current. For instance, to measure the piezoelectric effect, the sample is 
typically strained along the piezoelectric axis while being shorted out with 
a capacitor (see Fig. H). 

The theory discussed here only concerns phenomena where the macroscopic 
polarization is induced by a source other than an electric field. Even in this 
case, the polarization may (or may not) be accompanied by a field, depending 
on the boundary conditions chosen for the macroscopic sample. The theory 
addresses polarization differences in zero field: this concerns therefore lattice 
dynamics, piezoelectricity (as in the ideal experiment sketched in Fig. |^), and 
ferroelectricity. Notably, the theory reported here does not address the problem 
of evaluating the dielectric constant: this can be done using alternative 
approaches, such as the well-established linear-response theory |13|, |2^, or 



other more innovative theories 2S 



The bulk quantity of interest, to be compared with experimental measure- 
ments, is the polarization difference between two states of the given solid, 
connected by an adiabiatic transformation of the Hamiltonian. The electronic 
term in this difference is: 

/•At 

AP= dtJ{t), (14) 
Jo 

where J{t) is the current flowing through the sample while the potential is 
adiabatically varied, i.e. precisely the quantity discussed in the previous 
Section, Eq. ([TB|). Notice that in the adiabatic limit At goes to infinity and 
J{t) goes to zero, while Eq. ([T^) yields a finite value, which only depends on 
the initial and final states. We may therefore write: 

P = - lim — Im In (^|e^^^|^), (15) 

L^oo 2,11 

where it is understood that Eq. (|l^) is to be used twice, with the final and 
with the initial ground states, in order to evaluate the quantity of interest 
AP. Notice that L — >■ cx) in Eq. (plSj) is a rather unconventional limit, since 
the exponential operator goes formally to the identity, but the size of the 
system and the number of electrons in the wavefunction increase with L. 
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6 The case of independent electrons 



We now specialize to an uncorrelated system of independent electrons, whose 
A^-electron wavefunction is a Slater determinant. As discussed above, in 
this case the expectation value (X), Eq. (^2|), is uniquely determined by the 
one-body density matrix. However, the formulation is simpler when expressing 
{X) and the resulting polarization P directly in terms of the orbitals. 

We restore explicit spin variables from now on. Suppose is even, and 
|\E') is a singlet. The Slater determinant has thus the form: 

1^) = ;^^I<^1^1<^2^2 • • • ¥^N/2'^N/2\, (16) 

where cpi are the single-particle orbitals. It is then expedient to define 

|^)=e'¥^|^): (17) 



even l^f) is indeed a Slater determinant, where each orbital ipi{x) of l^') is 
multiplied by the plane wave e*~^. According to a well known theorem, 
the overlap amongst two determinants is equal to the determinant of the 
overlap matrix amongst the orbitals. We therefore define the matrix (of size 
N/2 X N/2y. 



-L 

Sij = {ipi\e'^''\ipj) = I dx^*(x)e^^>j(x), (18) 



in terms of which we easily get 

P = — ^Im In (^1^) = --Im In det S, (19) 

where the factor of 2 accounts for double spin occupancy, and the expression 
becomes accurate in the limit of a large system. 

The expression of Eq. (|19D goes under the name of "single-point Berry 
phase" (almost an oxymoron!), and was first proposed by the present author in 
a volume of lecture notes Since then, its three-dimensional generalization 



has been used in a series of DFT calculations for noncrystalline systems 
and has been scrutinized in some detail in Ref. |21[ . 



The case of a crystalline system of independent electrons is the one which 
historically has been solved first @, §], though along a very different logical 
path [0, 0, than adopted here. I am going to outline how the present 



formalism leads to the earlier results. 

For the sake of simplifying notations, I am going to consider the case of 
an insulator having only one completely occupied band. The single-particle 
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orbitals may be chosen in the Bloch form: indeed, the canonical ones must 
have the Bloch form. The many-body wavefunction, Eq. (0), becomes in the 
crystalline case: 



1^) = -^\Mq^Mq2 ■ ■ ■ ^giV/2^g^/J- (20) 



The lattice constant is a = 2L/N, and the Bloch vectors entering Eq. (pO|) are 
equally spaced in the reciprocal cell (0,27r/a]: 

47r 

gs = ^s, s = l,2,...,iV/2. (21) 

Owing to the orthogonality properties of the Bloch functions, the overlap 
matrix elements in Eq. (|18|) vanish except when g^' = Qs — 27r/L, that is 
s' = s — 1: therefore the determinant in Eq. ([T9|) factors as 



^ N/2 

P = -Imln n(^'?>'^1^'^-)' (22) 



s=l 



where V'9o(^) = V''?jv/2(^) implicitly understood (so-called periodic gauge). 

The expression in Eq. ( P^ is precisely the one first proposed by King- 
Smith and Vanderbilt p[ as a discretized form for the Berry phase: in fact, 
the multi-band three-dimensional generalization of Eq. (|22|) is the standard 
formula 0, |T3[ implemented in first-principle calculations of macroscopic 
polarization in crystalline dielectrics, either within DFT or HP. 



7 Critical rethinking of DFT 

The modern viewpoint about macroscopic polarization has even spawned a 
critical rethinking of density-functional theory in extended systems. The 
debate started in 1995 with a paper by Gonze, Ghosez, and Godby |2^, and 



continues these days |3y]. The treatment of polarization provided in the present 
work allows discussing the issue in a very simple way. 

The celebrated Hohenberg-Kohn (HK) theorem, upon which DFT is 



founded states that there exists a universal functional F[n], which de- 
termines the exact ground-state energy and all other ground-state properties 
of the system. The main hypotheses are that the magnetic field is vanish- 
ing, the ground state is non degenerate, and — most important to the present 
purposes — the system is finite, with a square-integrable ground wavefunction. 
This is precisely the key point when dealing with an extended system. 
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Ideally, it is possible to refer to a macroscopic but finite system: polariza- 
tion is then by definition the dipole divided by the volume, and Eq. (|I]) safely 
applies. As a consequence, the polarization of the real interacting system is 
identical to the one of the fictitious noninteracting Kohn-Sham (KS) system. 
Unfortunately, such polarization depends on the charge distribution both in 
the bulk and at the surface: according to Ref. ||2^, this fact implies a possible 
"ultranonlocality" in the KS potential: two systems having the same density 
in the bulk region may have qualitatively different KS potentials in the same 
region, and different polarizations as well. The issue can be formulated in a 
more transparent way by recasting it in the language of the present work. 

As discussed in the Introduction, condensed matter theory invariably works 
in a different way: one adopts BvK since the very beginning, and the system 
has no surface by construction. The original HK theorem was formulated 
for the case where the Schrodinger's equation is solved imposing square- 
integrable boundary conditions, but the same theorem holds within BvK, with 
an identical proof, for a finite A^-electron system. We can then define, even 
within BvK, the fictitious KS system of noninteracting electrons, having the 
same density as the interacting one: if the system is crystalline, the KS orbitals 
have the Bloch form. 

Now the question becomes: do the interacting system and the correspond- 
ing KS noninteracting one have the same polarization? The answer is actually 
"no": in fact, within BvK, polarization is not a function of the density, not 
even of the one-body density matrix, as stressed above. Numerical evidence 



of the fact that the two polarizations are not equal has been given |^ . 

Other important implications concern the occurrence of macroscopic 
electric fields within DFT: we refer to the original literature p9| , pO[| about this 
issue, while here we limit ourselves to just remarking an important point. The 
potential (both one-body and two-body) within the Schrodinger's equation 
must be BvK periodic, otherwise the Hamiltonian is an ill-defined operator in 
the Hilbert space. The periodicity of the potential is tantamount to enforcing 
a vanishing macroscopic electric field: therefore some ad-hoc strategies must 
be devised in order to cope with nonzero electric fields, as is indeed done in 



Refs. m, 30 
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